Full Spin and Spatial Symmetry Adapted Technique for Correlated Electronic 
Hamiltonians: Application to an Icosahedral Cluster 



Shaon Sahoo ' Q and S. Ramasesha*' | 

a. Department of Physics, Indian Institute of Science, Bangalore 560012, India. 

b. Solid State & Structural Chemistry Unit, Indian Institute of Science, Bangalore 560012, India. 



Abstract 

One of the long standing problems in quantum chemistry had been the inability to exploit 
full spatial and spin symmetry of an electronic Hamiltonian belonging to a non-Abelian 
point group. Here we present a general technique which can utilize all the symmetries 
of an electronic (magnetic) Hamiltonian to obtain its full eigenvalue spectrum. This is a 
hybrid method based on Valence Bond basis and the basis of constant z-component of the 
total spin. This technique is applicable to systems with any point group symmetry and is 
easy to implement on a computer. We illustrate the power of the method by applying it to 
a model icosahedral half-filled electronic system. This model spans a huge Hilbert space 
(dimension 1,778,966) and in the largest non-Abelian point group. The C^o molecule has 
this symmetry and hence our calculation throw light on the higher energy excited states of 
the bucky ball. This method can also be utilized to study finite temperature properties of 
strongly correlated systems within an exact diagonalization approach. 
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1 Introduction 



One of the major goals of the electronic structure theory of molecules is the determination of the 
excited states and their properties. For studying the linear and nonlinear optical properties of a 
system, we need to obtain excited states of desired symmetries, while we need the full excitation 
spectrum to study finite temperature properties. A brute force diagonalization of the full system 
Hamiltonian is not feasible even for a moderately sized systems and even if we succeed in 
obtaining all the eigenstates, it is difficult to identify them with irreducible representations to 
which they belong. In cases where we can manage to obtain the low-lying eigenstates of the 
Hamiltonian, we may miss the states important for the desired purpose, since in correlated 
many-body Hamiltonians, there can be an unpredictable number of 'intruder' states between 
the ground and desired excited state. Utilizing the full spacial and spin symmetry (conservation 
of total spin and z-component of total spin) allows one to obtain several low-lying eigenstates 
in each spatial symmetry subspace for every total spin value and for many low-temperature 
static properties of a system, this will suffice. For the study of dynamic properties as well as 
finite temperature properties, we need to know the full eigen spectrum. Obtaining the full eigen 
spectrum for a large molecular system is, however, not feasible by any method at the present 
time. But, utilization of all the symmetries of a Hamiltonian allows extending dynamic and 
finite temperature properties to a slightly larger systems than what is feasible in the absence of 
a symmetry. 

Most electronic structure calculations start with molecular orbitals and account for cor- 
relation by employing a configuration interaction (CI) approach either in a perturbative or a 
variational scheme. However, even a restricted CI approach, involving only frontier orbitals, 
becomes too difficult to handle for large molecules [2J. We can circumvent this difficulty by re- 
sorting to model Hamiltonian. In some molecular systems, it is possible to identify a subsystem 
to which the important electronic excitations are confined. In such a situation, it is both advanta- 
geous and insightful to deal with model electronic Hamiltonians which describe the excitations 
in the subsystem. One such molecular system is the conjugated n system. 

The model Hamiltonian for describing conjugated system was first introduced by Hiickel 
and has mainly served pedagogical purpose in understanding the chemistry of conjugated sys- 
tems [3]. More realistic models which take into account electronic repulsions within the n 
system was introduced by Pariser and Parr 01 as well as by Pople independently in 1953. 
This and related models such as the Hubbard model [6|| have dominated the study of correlated 
electronic systems in chemistry and physics for almost half-a-century. These models consist of 
a one-electron Hamiltonian defined in the basis of site orbitals and whose matrix elements are 
non-vanishing along the diagonal as well as between orbitals on chemically bonded sites and 
a two electron term which is approximated within a zero differential overlap (ZDO) scheme 
B71[8l. The ZDO scheme leads to electron repulsion integrals which are diagonal in the atomic 
orbital basis. The Pariser- Parr-Pople (PPP) model Hamiltonian is given by 
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Here, the first term of the Hamiltonian is the Hiickel term with c-^ (c i(J ) creating (annihilating) 
an electron of spin a at the i th site and the summation over bonded pair of sites < ij >. The 
second term is the Hubbard term with £/, being the on-site repulsion energy for i-th site (/?,- is 
the number operator for i site). The last part is the inter-site interaction term with Vy being 
the density-density electron-repulsion integral between sites i and j, Zj is the local chemical 
potential and corresponds to the occupancy of / site for which the site is neutral. We employ 
the Ohno interpolation scheme to parametrize Vy flU. 



Here ry is the distance (in A unit) between the i and /" sites using the Hubbard t/'s (in eV) 
at these sites. 

The Fock space of the PPP Hamiltonian scales as 4 N where N is the number of orbitals 
considered in the system and obtaining even a few exact low-lying states of the Hamiltonian for 
reasonable N could pose a challenge. While this problem can be managed to some extent by 
resorting to approximate treatments such as restricted CI schemes by (1) restricting the number 
of active orbitals considered in the CI step and (2) by considering only some classes of particle- 
hole excitations of the system [2], the advantage of exploiting all the symmetries possessed by 
the PPP Hamiltonian cannot be overstated. Full symmetry adaptation, besides factorizing the 
Hilbert space and thereby reducing computational effort also provides the symmetry labels of 
the states for discerning the state properties. The PPP Hamiltonian, being non-relativistic con- 
serves total spin, S, as well as z-component of total spin, Ms and could possess additional spatial 
symmetries depending on the system in question. The diagonalization of the Hamiltonian can 
be simplified by specializing the basis, in which the matrix representation of the Hamiltonian is 
sought, to the case of fixed total spin and z-component of the total spin and a specific irreducible 
representation of the point group. 

The conservation of the Sj ot , the total z-component of spin is achieved by choosing from the 
Fock space, states whose total Ms corresponds to the desired value. This is trivially possible by 
choosing a spin orbital basis and populating them with electrons to obtain the desired total Ms- 
It is also quite straightforward to set up the Hamiltonian matrix in this basis and solve for a few 
low-lying states in cases where the Hilbert space is spanned by a few hundred million states 
(see subsection 2.1). Factorizing the Hilbert space into different irreducible representations of 
the point group of the Hamiltonian is also straightforward as the resultant of a spatial symmetry 
operator, operating on a Slater determinant is easy to obtain in atomic orbital basis. In modern 
quantum chemical calculations, these symmetries are routinely employed. 
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However, construction of spin adapted configuration state functions which are simultaneous 
eigenstates of Sf ot and Sj ot operators is nontrivial and pursuit of this has been a long standing 
interest in quantum chemistry. The Hamiltonian matrix in such a symmetrized basis leads to 
matrices of smaller order besides allowing automatic labeling of the states by the total spin. 
Furthermore, the eigenvalue spectrum is enriched, since we can obtain several low-lying states 
in each total spin sector. This can be contrasted with obtaining several low-lying states in a 
given total Ms sector which would have states with total spin S tot > Ms- There are many ways 
of achieving this task ifTOl : most important among these are Valence Bond (VB) approach ifTTI . 
Lowdin spin projection technique lfT2l [T3l and group theoretical approaches lfT4l [T5l . While 
they are satisfactory regarding spin adaptation, most of these techniques virtually fail while 
dealing with non-Abelian spatial symmetry. They become symmetry-specific, even frequently 
impractical while applied to large system with a non-Abelian symmetry (see review in [Q]|). 
Here we present our hybrid VB-constant Ms method, which overcomes these difficulties. 

The ultimate goal of symmetry adaptation is to exploit the full spatial and spin symmetries of 
the system, both for computational efficiency and for complete labeling of an eigenstate by the 
total spin and the irreducible representation it which belongs. In Sec. 2, we present our hybrid 
VB-constant Ms method which allows exploiting the full spin and spatial symmetries of any 
arbitrary point group. Similar method applicable only to pure spin systems has recently been 
developed [1J. The technique presented here is applicable to more general systems of correlated 
electrons. In Sec. 3, we illustrate an application of this method to a PPP and Hubbard model 
of the half-filled icosahedron which has one orbital at each of the 12 vertexes. The icosahedron 
is the smallest system with all the symmetries of C^q, the carbon Bucky ball and obtaining all 
the eigenstates of this model will throw light on the correlated states of C(,q. In Sec. 4, we 
summarize and discuss the technique. 

2 Hybrid VB and Constant M s Basis Method 

In an electronic system, a given orbital can be in one of four states; it can be (i) empty, (ii)singly 
occupied with an up spin electron, (iii) singly occupied with a down spin electron and (iv) 
can be doubly occupied. Constant Ms bases, for a given filling of the orbitals, are obtained 
trivially by choosing states from Fock space, whose total S z value corresponds to the desired 
Ms value (Ms = sum of z-components of individual electron-spins). By construction they are 
orthonormal. The easiest way of constructing the spin adapted functions is the diagrammatic 
valence bond (VB) method based on Rumer-Pauling rules |fT6l ITTTl . If N is the number of 
orbitals, N e is the number of electrons with up-spin electrons and N± down-spin electrons 
(N e = +N\), then, all possible linearly independent and complete set of states with total spin 
S and Ms = S, for a fixed occupancy of the orbitals, according to extended Rumer-Pauling rules 
are obtained as follows, (i) The N orbitals are arranged as dots on a straight line, (i) Doubly 
occupied sites are marked as crosses, (ii) An arrow is passed through 2S of the singly occupied 
vertexes, passing on or above the straight line on which the system is represented. The arrow 
denotes the spin coupling corresponding to total spin S and total z-component Ms = S. (iii) 
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Figure 1 : Representation of VB diagrams for a half-filled 6 orbital (site) system. Here • denotes empty 
site, x denotes doubly occupied site. The top VB diagram shows a spin pairings to yield a state with 
total spin 5 fot =l, its bit representation corresponds to a unique integer I = 2350. The bottom VB diagram 
shows a Stat = state, the corresponding unique integer, I is 2169. 



Remaining singly occupied vertexes are singlet paired and are denoted by lines drawn between 
them which lie on or above the straight line describing the system, (iv) Diagrams with (a) two 
or more crossing lines or (b) crossing line and the arrow or (c) a line enclosing the arrow are 
rejected. The remaining set of diagrams correspond to a complete and linearly independent 
set of VB states for the chosen orbital occupancy. The set of VB diagrams which obey the 
extended Rumer-Pauling rules would hence forth be called "legal" VB diagrams. Some legal 
VB diagrams are shown in Fig. (OQ) along with the integers which represent them. In the case 
of N e odd and S = 1/2, we cannot have an arrow with just one site! We handle this situation 
by augmenting the system by adding a "phantom" site. The VB states of all legal singlets with 
single occupancy of the phantom site provides the complete and linearly independent basis. The 
phantom site appears only in the basis and not in the system Hamiltonian. 

We can generate the complete set of VB states for our case of N orbitals with N e electrons 
of total spin S and z-component of total spin Ms = S by exhausting all possible occupancies of 
orbitals which satisfy 2S = N+ — N\. Since an orbital can be in any one of four states (empty, 
doubly occupied, a singlet line beginning or a singlet line ending, sites involved in an arrow 
being treated as line beginnings) we can use two bits to represent the state of an orbital. Thus, 
each VB diagram can be uniquely represented as an integer on a computer. 

A line in the VB diagram, between sites "/" and "/', represents {a] a a\ n — fijadj a ) |0 >/\/2, 
where we choose a to correspond to | f) and |3 to 1 1) orientations of the electron. The doubly 
occupied site "/" corresponds to the state aj a dj^\0 >. The phase convention assumed for a line 
between sites "z" and "/' is that the ordinal number "f ' is less than the ordinal number "/'. 
The 25 singly occupied sites k\ £2 £3 . . . . &2S m the arrow represent the state with Ms = S 
given by a aJ 2 a aJ 3 a ---«J 25 a |0 >. VB states corresponding to other M5 value for this state 
with spin S, can be obtained by operating, required number of times by the S wt operator on the 
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state. Since Hamiltonian in Eq.® is isotropic, each eigenstate in the spin S sector is (2S+1) fold 
degenerate and, by Wigner-Eckart theorem IfTTll it is sufficient to work in subspace of chosen 
Ms value. The VB state corresponding to a given diagram is a product of the states representing 
the constituent parts of the diagram, in no particular order as each part is either a product of two 
Fermion operators or a linear combination of the product of two Fermion operators. 

Given the definition of a line in the VB diagram, every VB diagram, |\]/,- > can be broken up 
into a linear combination of the constant Ms basis states >} as, 

|Yi>=£Cy|4>y> (3) 

j 

A VB diagram with n singlet lines yields 2" basis states in the constant Ms basis. To effect 
the conversion of VB diagrams to constant Ms functions, we note that each singlet line gives 
two states; in one state, the site at which a singlet line begins is replaced by an a spin while the 
one at which it ends by a P spin with phase +1 and in the other the spins are reversed and the 
associated phase is -1. There is a normalization constant, (2~"/ 2 , associated with the constant 
Ms basis state. The matrix relating the VB basis states to constant Ms basis states, C, is a V x M 
matrix, where V is the dimensionality of the VB space and M that of the constant Ms space. If 
Rm the matrix representation of symmetry operation R is known, in constant Ms basis, then the 
knowledge of the matrices C and gives the result of operating by the symmetry operator 
R on a VB state as a linear combination of the constant Ms basis states via the matrix = 
CRm- The projection operator for projecting out the basis states on to a chosen irreducible 
representation T of the point group is given by, 

Pr = X>F(*)tf (4) 

R 

where, %f r (R) 1S me character under the symmetry operation R of the point group of the system 
fl2TTl . The matrix representation of Pr in the mixed VB and constant Ms basis is given by, 

Qr = £xf r (^)B J! (5) 

R 

where, Qr is a V x M matrix. However, the rows of the matrix Qr are not linearly independent, 
since the complete symmetrized basis transforming as T spans a much smaller dimensional 
Hilbert space. The exact dimension Vr of the Hilbert space spanned by the system in the 
irreducible representation T can be known a priori and is given by, 

V r = (d r /h)£ X red (R)Xr r (R) (6) 

R 
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where dr is the dimensionality of the irreducible representation T, h is the number of symmetry 
elements in the point group and % red (R) is the reducible character for the operation R. The 
determination of % red (R) is nontrivial and the method of computing it will be discussed in the 
next subsection. The Vr x M projection matrix, Pp of rank Vr is obtained by Gramm-Schmidt 
orthonormalization of the rows of the matrix Qp until Vp orthonormal rows are obtained. These 
orthonormal and linearly independent rows yield the desired linear combinations which trans- 
form as r and also have total spin S. Projection matrix Pr is represented by these Vr rows. 

The M x M Hamiltonian matrix is constructed in the constant Ms basis (see subsection 
2.1). Since the basis states in this representation are orthonormal, we do not encounter the 
problem of "illegal" VB states. In the pure VB method, the Hamiltonian operating on a legal 
VB state can yield illegal VB diagrams which then need to be re-expressed as linear combination 
of the legal VB functions [fT9l . The Vr x Vr Hamiltonian matrix in the fully symmetrized basis 
is given by PrHw P r and one could use any of the well known full diagonalization routines to 
obtain the full eigenspectrum or use the Rettrup modification of Davidson algorithm [|22l to get 
a few low-lying states of the symmetrized block Hamiltonian in the chosen spin and symmetry 
subspace. 

For degenerate irreducible representations, such as the E, T, G or H representations, the 
above procedure does not lead to the smallest block of the Hamiltonian matrix. In such cases, 
it is advantageous to work with bases that transform according to one of the components of the 
irreducible representation. In case of E, T or H, this can be achieved by choosing an axis of 
quantization and projecting out basis states of the irreducible representation which are diagonal 
about a rotation about the quantization axes. For example, in the case of the irreducible rep- 
resentation that transforms as T, we can choose one of the C3 axes as a quantization axis and 
project the basis states which transform as the irreducible representation T, using (I + C\ + C|) 
as the projection operator. This operator projects states that transform as the Y® component of 
the three fold degenerate irreducible tensor operator. Similarly we can choose a C5 axis and use 
(7 + C\ + Cg + Cg + Cg) as projection operator for the irreducible representation H. For the E 
representation, we can use (I + C2) as projection operator with a chosen C 2 axis. The case of G 
is a bit tricky; one has to choose two C 2 axes, orthogonal to each other. The projection operator 
for this case then would be: (I + C2) (I + C' 2 ) ■ After these projection operations, dimensions of 
the Hamiltonians to be diagonalized would be half, one third, one fourth or one fifth respectively 
for the E, T, G and H representations. 

Here we wish to emphasize the computational advantage of our technique over the constant 
Ms basis method. The additional steps involved in the hybrid VB-Constant M5 method are 
(i) construction of the C matrix and (ii) computation of the matrix. However, if we wish to 
compute the properties of a state expressed as a linear combination of VB diagrams, the simplest 
way is to use the C matrix to transform the state from the VB basis to the constant Ms basis. 
Therefore, construction of the C matrix is not strictly an overhead. Besides,the construction of 
the C matrix is a very fast step as the row index of the C matrix is the index of the VB state 
which we wish to decompose and the column indices of elements in this row are the indices 
of the resultant constant M5 states. The constant M5 states are easily generated as an ordered 
sequence of integers which represent them and this facilitates searching for the column index of 
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the matrix. The coefficients will have a magnitude of 2~ n l 2 where n is the number of lines in the 
jth yg di a g ram ; the phase of the coefficient is easily fixed based on the phase convention used 
for a singlet line. In the hybrid approach, computation of the matrix involves the matrix 
multiplication, CRm- The number of arithmetic operations involved is however very small, 
since both C and are sparse matrices with the latter having only one nonzero matrix element 
per row. In both constant M5 and hybrid approaches one has to obtain the projection matrix Pp 
by retaining only the orthogonal rows of the matrix Qr- Since the number of orthogonal rows 
in Qr is far fewer than in Rm, this step is faster in the hybrid approach than in the constant 
Ms approach by a factor D(Ts)/D(Tm s ), where D(Ts) is the dimensionality of the space of the 
irreducible representation T with spin S and D(Tm s ) is similarly the dimension of the space T 
with constant M5. Though, this advantage is largely off-set by the fact that the R^ matrix in 
constant Ms basis is more sparse than the Qr matrix in the hybrid approach. Computation of 
the eigenvalues (diagonalization of Hamiltonian) in the constant Ms approach is slower than in 
the hybrid approach, since D(rM s )>D(T5) for most S (for example see TableQ])- The memory 
required for the hybrid approach is not very different from that of constant Ms approach, even 
though the matrices in the hybrid approach are slightly denser, they are smaller in size. The 
only additional memory demand in the hybrid approach is the storage of sparse C matrix. The 
major advantage of the hybrid approach is that we can obtain a far richer spectrum, since we 
are targeting each spin sector separately, unlike in the constant Ms approach. Thus, if we can 
obtain (by our approach), say 10 states in each S sector, then each one will correspond to a 
unique state. There will be no repetition of the states. But in contrast, 10 states obtained in a 
Ms sector (by constant Ms approach) may not be unique, since many of these states would be 
repeated in different Ms sectors. 

2.1 Implementation Details 

We can represent a basis uniquely by a 2N-digit binary number [11] (N is the number of sites 
/ orbitals); the first two bits describe the state of the first site, next two bits describe the state 
of the second and so on. For constant Ms basis, we use the bit states "00" for an empty site, 
"10" for site with spin-up electron, "01" for site with spin-down electron and "11" for a doubly 
occupied site. Similarly, for VB basis, we use "00" for empty site, "10" for a singlet line- 
beginning at a site as well as for all sites in the arrow, "01" for a line-ending and "11" for a 
doubly occupied site. The Rumer-Pauling rules are implemented by enforcing that for a given 
site the quantity (# of line beginnings - # of line endings) at sites one to i — 1 should be > 
|fT8l . Our binary coding implies that the number of bits in the state "1" is equal to the number 
of electrons N e . To generate integers that represent VB states, we generate integers with N e "1" 
bits in the bit field from zero to (2N-1) in increasing sequence and check to see if the bit pattern 
corresponds to a "legal" VB diagram with chosen total spin. For integers corresponding to 
constant Ms basis, total Ms value should be the desired value and the Rumer-Pauling condition 
is not enforced. The positive integers so generated uniquely represent the states of the VB or 
constant Ms bases. 

Computationally, finding the transformation matrix C which carries the VB basis to constant 
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Ms basis is straightforward. We initialize the coefficients in the row of the matrix C correspond- 
ing to the chosen VB state to zero. We then decompose the VB diagram by converting every 
singlet line in the diagram into two Ms states. The indices of the resulting constant Ms states 
correspond to the column indices of C and are determined by a binary search on the list of 
integers that represent the constant Ms states. The corresponding matrix element is given by the 
normalized VB coefficient with appropriate phase. On a computer, the transformation matrix 
C is stored in sparse form. Next, we construct the projection matrix (Pr), by constructing the 
matrix representation of each of the symmetry operators, R, of the point group in the constant 
Ms basis. This is achieved by (i) obtaining the occupancies of each site from the integer repre- 
senting the basis state, and (ii) by letting R act on the basis state by appropriately rearranging 
the sites together with their occupancies to obtain the new bit pattern corresponding to the re- 
sultant state. The new occupancy pattern is converted into the integer representing the state and 
fixing the column index of the matrix R by a binary search for the index of the new integer 
in the list of integers representing the constant Ms basis. Care should be taken to keep track 
of the phase factor while interchanging the occupancies since Fermion creation operators anti- 
commute. From a knowledge of C and all the Rm matrices we can construct the matrix Qr (Eq. 
[5]). But the rows of Qr are in general not linearly independent, eliminating linear dependencies 
leads to the projection matrix Pp with Vp number of linearly independent rows. The Vr linearly 
independent rows can be obtained by (i) collecting all linearly independent rows, by inspection, 
by noting that the set of rows which are disjoint (that is do not have non zero elements with 
common column index) are orthogonal by virtue of the fact that the constant Ms basis sets are 
orthogonal and (ii) by carrying out Gram-Schmidt orthonormalization to obtain the remaining 
linearly independent rows. However, knowing Vr a priori is important to be able to stop the 
orthonormalization process once the number of linearly independent rows obtained equals the 
dimensionality of the symmetrized space. While Vp can be obtained from Eq. [6l it needs a 
knowledge of the reducible character. 

To obtain the reducible character, it appears that we need a matrix representation of the 
symmetry operator in the VB basis. Given an operator R, the matrix representation r, in the VB 
basis, is obtained from, 



%/>=£ry|\|fy>. (7) 

j 

However, since the VB basis is non-orthogonal, we need the inverse of the overlap matrix,S _1 
where Sy =< > are tne matrix elements of S. The matrix r is then given by RS _1 , where 
the matrix elements of R are given by Rij =< \\fj\R\\\fj >. In general determination of the ma- 
trix S _1 is difficult for Fermionic systems and computationally prohibitive for large pure spin 
systems. 

The above difficulty can be circumvented by resorting to the bit representation of VB and 
constant Ms basis states. Using the C matrix, we can rewrite [7J as 
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j k 



(8) 



For every state >, we need to find the coefficient r ; , and the reducible character % red (R) = 
Y.i Hi- Taking the inner product on both sides of Eq. [8] with |(|)/ >,we get, 

< §i\R\Wi >= J^rijC Jk < §i\$k > (9) 

j k 
j 

rij are unknowns and need to be determined. The Ihs can also be evaluated as 

<Q> l \&\y i >=<$ l \R\ y £c ij \Q> j > (10) 

j 

j k 

= Y. R lkCil 

k 

where Rij is the matrix representation of R in the constant Ms basis which is known. The only 
unknowns on the rhs of Eq. [9] are the coefficients ry and we need to determine the diagonal 
elements r#. 

To determine r„, let us first assume that integers {7/} represent the constant Ms basis states 
>} and the integers {/,•} represent the VB states {|\|/; >}. Now we note that in the expansion 
of a VB state |\|/ ; > as a linear combination of constant Ms states (Eq. [3]), the largest integer that 
represents the constant Ms state, 7/ which appears in the expansion, is the one corresponding to 
the integer that represents the VB state > itself, namely This is because, we have chosen 
the bit state "10" both for a line beginning in the VB state and for an up spin occupancy in 
the constant Ms basis. The assertion that U > 7 7 in the decomposition of the VB state |\|/; > 
into constant Ms functions \tyj >, implies that the matrix C has nonzero elements Cjk only for 
Jk<Ij- 

Let us consider Eq. [71 we note that on the rhs the summation runs over all the states of the 
VB basis. Let us consider the VB state, |\|/y > which is represented by the largest permitted 
integer, Iy. This integer also correspond to the Ms basis state \§m > (where M is the dimen- 
sionality of the constant Ms space) with the largest integer representation, {Iy — 7m), Taking the 
inner product with the state \§m >> from Eq. [10] and Eq. (9[ we obtain, 



< § M \R\yi >= nyC v ,M 
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All other terms on the rhs of Eq. [9] are zero. Hence using Eq. QT1 we can determine r^y. 
We can now proceed with the constant Ms state whose representing integer J% is equal to Iy-i- 
The constant Ms state with Jk can appear only in the expansions of the VB states \|/y and \\fy-i . 
Taking the inner product with we obtain, 

< §K\R\Wi >= nyC V ,K + n,V-lC V -LK (12) 

In [T2l the only unknown is r^y-i and can be evaluated. Similarly, by proceeding to the 
VB state \|/y-2, we can obtain r,y_2- We can terminate when we reach the VB state |\|/; >. 
This procedure can be adopted to obtain all the diagonal elements of the r matrix and hence the 
reducible character, % red . 

Constructing the Hamiltonian matrix in constant Ms basis in real space is a fast and easy 
step. The basis states are eigenstates of the interaction part for the model Hamiltonians. From 
the binary sequence of the integers which represent the constant Ms basis, we know the occu- 
pancy of each site and hence the diagonal contribution of the interaction terms. Simple rules 
for operating on a constant Ms basis state by the operators Eij = £ (aJ a /0 + at o a (a ) have been 
published elsewhere and together with a binary search procedure which allows rapid generation 
of the matrix corresponding to the one-electron terms of the Hamiltonian [[TT|. 

3 Application to Icosahedral Cluster 

To illustrate the power of our technique, we have applied the method to a 12-site regular icosa- 
hedral cluster (see Fig. [2]) at half-filling. It has 30 edges (each one here taken to be of length 1.4 
A) representing a chemical bond.We have chosen this system, because it belongs to very high 
symmetry non-Abelian point group and presents a very general case for testing our method. 
This point group is also the same as the point group of the C(,o molecule. So properties, which 
are particularly symmetry-related, obtained for our model system, would also be useful in gain- 
ing insights into the C^q molecule. 

We have studied the icosahedron within the Hubbard model in which the inter-site interac- 
tions are neglected, for a range of U values, as well as in the PPP model with standard Carbon 
parameters. The number of Ms=0 states is 853,776 and we have obtained the exact energies 
of all the states, by using the full spatial and spin symmetries of the system. The energies of 
Ms 7^ states are also known, since we know the total spin of each state. We have studied the 
density of states in each symmetry and spin sector as a function of U in the Hubbard models 
and also for standard parameters in the PPP model. 

In Table CQ), we give the dimensions of all the subspaces of different total spin and total Ms 
values, for the Icosahedral cluster. 

We note here the huge fall in size of total spin subspaces compared to total Ms subspaces 
for most of the cases. Using the hybrid VB-constant Ms method, we have broken down each 
total spin sectors into basis states that transform as different irreducible representations of the 
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Figure 2: Regular icosahedron with (//,) symmetty. Our system has an orbital at each of the 12 vertexes, 
and a transfer integral corresponding to a bond on each of the 30 edges; In the PPP model, the bond 
length is taken to be 1.4 A. 



icosahedral point group. The dimensionalities of the various symmetry subspaces are shown in 
Table©. 

We note that most subspaces are small enough for obtaining all the eigenstates of the Hamil- 
tonian. However, for degenerate representations the subspaces are large and can be reduced by 
a factor equal to the dimensionality of the representation, as described earlier. We have used 
this approach and obtained all the eigenstates of the ensuing Hamiltonian matrix using a full 
matrix diagonalization routine. Since the number of eigenstates in each subspace is large, we 
have computed the density of states (DoS) using a AE of 0.4eV, for which the histograms of the 
DoS are stable. A histogram for particular spin is evaluated by summing over corresponding 

Table 1 : Dimensionalities of different spin subspaces of a half-filled 12-site (orbital) electronic system. 
D(S) is the dimensionality of the constant S basis and D(Ms) is the dimensionality of the constant Ms 
basis. 
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Table 2: Dimensionalities of different symmetry and spin subspaces of half-filled icosahedral cluster. 
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50932 
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16 





Hg 


47940 
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9220 


900 


40 





K 
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3188 


1644 
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40 








Ti„ 


17082 


28686 


14700 


3372 


294 
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17082 


28686 


14700 


3372 


294 


18 







30160 


50992 


26176 


5888 


560 


16 





H u 


46880 


79680 


40980 


9060 


900 


20 





Tot Dim -> 


226512 


382239 


196625 


44044 


4212 


143 


1 



states of all symmetries. Same is for histogram for particular symmetry, where corresponding 
states of all spins are considered. Although, unlike the one-particle DoS, the many-body DoS 
is not an intensive quantity. However, for a given model and system size, we can use these 
quantities to understand the behavior of the system. 

3.1 Hubbard model studies 

In Fig. ©, we show the many-body DoS for the Hiickel model in the the various A g , H g and 
A u , H u spaces. Here, we have summed over all spin states, for simplicity. We have also shown 
the DoS plots of different spin space, in which they are summed over all the irreducible spaces. 
Two things are worth noting. Firstly, the DoS displays a symmetry about zero of energy in 
each of the subspaces, even though the system does not possess the e-h symmetry. In fact, it 
is clear from the one particle spectrum that there is no symmetry in the one-particle energy 
levels about zero energy. Secondly, the DoS profile of the g subspace shows peaks where 
ever there is a valley in the DoS profile of the u subspace. This is indeed true also for other 
irreducible representations not shown in the figure. The reason for this symmetry in the DoS 
plots is because the sum of the one particle eigenvalues are zero and follows from the fact that 
in Hiickel model, with all site energies set to zero, the diagonal matrix elements are all zero. 
This implies £,2e, = 0, where £, are the molecular orbitals (MO) energies. Thus, for any 
given occupancy pattern of the MOs at half-filling, we have J^n,-£j = — £, (2 — «,-)£,- and since 
L; n i — Hi (2 — ni), at half filling, we find that for every many-body state of energy there exists 
a many-body state of energy — E^, even though the MO energies £/ do not satisfy the pairing 
theorem [20J. Thus, the symmetry in DoS plots is not a consequence of the pairing theorem 
but due to the fact that the magnitude of the sum of the energies of the bonding MOs is equal 
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Figure 3: DoS profiles for (a) A g and A u and (b) H g and H l( spaces. For other symmetry subspaces, the 
DoS profiles are similar. We note the difference in profiles for g and u spaces. DoS profiles for (c) S = 
and (d) S = 1 are also given. Inset of (b) gives the one-particle spectrum for the regular icosahedron in 
the Hiickel model. 

to the magnitude of the sum of the energies of anti-bonding MOs. This is in fact a general 
result for the Hiickel model with equivalent sites. The second observation about the location 
of valleys and peaks in the DoS of the states with g and u symmetries is due to the fact that 
the molecular orbital occupancies which give the g and u representations are different due to 
symmetry considerations. The DoS plots in various spin subspaces are also shown in Fig. ©. 
We find that they show several peaks in each total spin space. We show in Fig. © DoS profiles 
for the Hubbard model in various symmetry subspaces for different values of \U Jt\ and in Fig. 
([5]) we show the same in different total spin spaces. The evolution of DoS with correlation 
strength is interesting. Firstly, we note that for \ U/t\ = 2.0, the sharp peaks in the DoS found in 
the Hiickel model are broadened. The peaks in the g subspace coincide with the troughs in the 
u symmetry and vice versa. The ground state energy in the presence of correlations is higher 
than in the Hiickel model, as expected. 

In the very strong correlation regime, (\U/t\ = 12 and \U/t\ =40 we again find peaks in 
the DoS in all the subspaces. What is interesting is that peaks appear at almost the same value 
of energy in g and u subspaces, unlike in the non-interacting or weakly interacting model, and 
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Figure 4: DoS profiles in various symmetry subspaces for four different \U ft\ values are shown. Note, 
in (a), profiles for T\ g and T2g coincide and simply referred to as T g . Same is for the triply degenerate u 
space. In (b), (c) and (d), DoS profiles for corresponding u and g spaces coincide, so they are referred to 
by their common irreducible representation symbols. 



are approximately \U/t\ apart in energy. This can be understood by noting that the many body 
space can be subdivided into space of all singly occupied sites; space of one empty, one doubly 
occupied and rest singly occupied sites; space of two empty, two doubly occupied and rest 
singly occupied sites and so on. The interaction energies of these class of states is 0, U, 2U, 
etc. The transfer term leads to weak admixture of these states, and in the strong correlation 
limit results in broadening of the DoS peaks centered at the energies 0, U, 2U etc. Thus, the 
DoS plots, although look similar in both the small \U/t\ and the large \U/t\ limits, their origin 
as well as their location is different. It is also worth noting that the DoS curves centered at 
different energies are similar for all the subspaces. However, in the large \U/t\ limit, the DoS 
curves for the same total spin centered around different energies are not similar, showing that 
we do not have a strict spin-charge separation in icosahedron in this limit. For, if indeed we 
had such a separation, we would expect very similar DoS for the same total spin, for different 
number of doubly occupied sites, when the allowed number of total spin states is large. 

For \U/t\ = 8.0, the DoS is very different. We find that there is a single broad peak and all 
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Figure 5: DoS profiles of various spin subspaces for four different \U /t\ values are shown. 



the other peaks are small inflexions superimposed over the peak. The "band width" of the one- 
particle spectrum (see Fig. [3b, inset) is 7.236? and the Hubbard correlation strength is very close 
to this value. Thus, for an icosahedral cluster, the parameters are at the intermediate correlation 
regime and leads to a smearing of the structure which is found at the weak and strong correlation 
limits. However, independent of the correlation strength, we find that the DoS curves are nearly 
identical for the Ti g and T2 g spaces and also for the T\ u and T% u spaces. 



3.2 PPP model studies 

The Hubbard model is not the appropriate model for studying carbon systems as it neglects 
long-range interactions. The appropriate model for studying such systems is the PPP model, 
which we have employed for studying the cluster. 

In Table ©, we show two lowest energy levels in each of the subspaces. We note that the 
ground state is the lowest energy state in the A g subspace with total spin zero. The one-photon 
gap is given by the lowest energy excitation to the T\ u space for an Icosahedron. Thus, we find 
that lowest excitation gap is at energy of 3.846 eV. This can be compared with the excitation gap 
of 3.552 eV for a PPP chain of 12 carbon atoms 11231 . We do not compare the excitation gaps to 
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Table 3: Lowest and second lowest energies for each symmetrized spin sector (for PPP model). All 
energies are in eV. 
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26.206 




0.831 


1.393 


5.401 
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PPP ring of 12 carbon atoms, as the N = An (n integer) have very strong finite size effects due 
to degenerate partly filled highest occupied molecular orbitals. The second allowed excitation 
is at an energy of 4.051 eV. The two photon gap is to the second lowest energy state in the A g 
representation and is found to be 1.533 eV which is very low compared to a polyene chain (~ 
3.0 eV). 

The lowest energy spin gap, from the singlet ground state is to the lowest energy spin 1 
state in the H g space. This gap is 0.122 eV which is very low compared to the polyenes. 
In general, the singlet-triplet gaps in conjugated systems is about 60% the optical gap and 
icosahedron seems to be an exception. Since singlet-triplet gap here is much higher compared 
to room temperature enery (about 0.025 eV), the system would show diamagnetic behavior 
below this temperature. There is also another triplet state which is about 0.012 eV above the 
lowest energy triplet state. These observations imply that the icosahedral cluster would exhibit 
paramagnetism above room temperature due to significant population of these states. Based on 
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Figure 6: DoS profiles of various symmetry subspaces for PPP model. 



the similar argument, we also conclude that the specific heat at low-temperature will be very 
small and increase exponentially with increasing temperature. The triplet-triplet (TT) excitation 
from the H g space is to states in T\ u , T2 U , G u and H u while from the G g state is to states T2 U , 
G u and H u . This means that we would have a band of TT excitations starting from 2.11 eV. 
Regarding higher spin excitations, there is a low energy quintet state about 0.279 eV above the 
ground state and a few other quintet excitations of energies between 2.414 and 2.937 eV. All 
other spin excitations are very high energy excitations, as the higher spin states have very low 
kinetic stabilization. 

In Figs. © and CO), we show the density of states plots for different symmetry subspaces 
and different total spins. We note from the figures that the icosahedral cluster of conjugated 
Carbon atoms belongs to the weakly correlated regime since the DoS peaks in the g and u 
spaces do not coincide in energy. We also find that this conclusion is corroborated by the DoS 
plots for different total spin states. The 5 = DoS plot shows a featureless broad peak, as seen 
for small \U/t\ values of the Hubbard model. The higher spin states also show broad peaks 
consistent with the weak correlation regime. Even in the PPP model, the DoS plots for the 7i g 
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Figure 7: DoS profiles of various spin subspaces for PPP model. 



and T2 g are very similar and same is the case with the T\ u and Ti u states. These DoS plots would 
also indicate the nature of the electronic spectra in these systems. 

One of the most fascinating molecules to have been discovered is C6o, which also has icosa- 
hedral symmetry. The Hiickel band width of C(,q is 5.618 t, when the transfer integrals for 
the hexagon-pentagon and the hexagon-hexagon bonds are taken to be the same. This is much 
smaller than the 7.236 t found for icosahedron. This is largely due to the different number of 
bonds per site (2.5 bonds / site for icosahedron compared to 1.5 bonds / site for C6o) in the 
two systems. For this reason, we expect PPP model with standard parameters of C6o to be in a 
more strongly correlated regime than the icosahedron. This should also reflect in the electronic 
spectra of Cgo- The full spectrum of icosahedron will also be helpful in gaining insights into the 
contributions of different states to the linear and nonlinear optical response of the system. 
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4 Summary 



In this paper, we have considered the long standing problem of both spacial and spin symmetry 
adaptation for arbitrary point groups. We have shown that by using the strengths of the VB 
and the constant Ms methods, we can have a hybrid scheme which exploits the full symmetry 
of a non-relativistic Hamiltonian. We have illustrated this by applying to the nontrivial case of 
an icosahedral cluster. We have obtained all the eigenstates of the cluster by our method. The 
hybrid method is less demanding on both memory and CPU time of a computer and is easy to 
implement. We have obtained the DoS of the Hubbard model for different Hubbard parameters 
U and of the PPP model for an icosahedral cluster. These plots show different characteristics 
as a function of interaction strength in the Hubbard model. The PPP model studies indicate 
that while the one-photon gaps are comparable with other conjugated systems, the spin gaps 
are unusually small. This may lead to significant population of the magnetic states at room 
temperature. These studies have a bearing on the system which also possesses icosahedral 
symmetry. The method discussed here will be of considerable importance in studying the dy- 
namics and finite temperature properties of systems whose Hamiltonians are amenable to exact 
diagonalization. While we have illustrated the method using a highly symmetric Hamiltonian, 
the method is very general and applicable to systems belonging to any point group. In point 
groups with lower symmetry, while the advantage of automatically labeling the states exists, 
the actual savings in computational effort would be decreased. 
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